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Abstract 

We compute a number of universal amplitude ratios in the three-dimensional Ising universal- 
ity class. To this end, we perform Monte Carlo simulations of the improved Blume-Capel model 
on the simple cubic lattice. For example, we obtain A^/A^ = 0.536(2) and C+/C- = 4.713(7), 
where A± and C± are the amplitudes of the specific heat and the magnetic susceptibility, respec- 
tively. The subscripts + and — indicate the high and the low temperature phase, respectively. 
We compare our results with those obtained from previous Monte Carlo simulations, high and 
low temperature series expansions, field theoretic methods and experiments. 
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I. INTRODUCTION 



In the neighborhood of a second order phase transition various quantities diverge, 
following power laws. For example, in a magnetic system, the correlation length ^, the 
magnetic susceptibility x ^'^^ the specific heat C behave as 

e^/±itr^ , x^c^w ,c^A±|tr-, (i) 

where t = {T — Tc)/Tc is the reduced temperature. The symbol ~ means asymptotically 
equal; corrections vanish as t — )■ 0. Critical exponents like u, 7 and a are universal. That 
means, they take exactly the same value for all systems in a given universality class. A 
universality class is characterized by the spacial dimension of the system, the range of the 
interaction and the symmetry of the order parameter. For reviews on critical phenomena 
and its modern theory, the renormalization group see for example 

While individual amplitudes like /+, /_, C+, C_, and A_ depend on the details of 
the system, amplitude ratios like f^/f^, C^/C^ and Aj^/A_ are universal. The indices + 
and — indicate the high and the low temperature phase, respectively. In addition to these 
simple ratios, there are also more complicated universal combinations of amplitudes. The 
combinations of the corresponding quantities are dimensionless. This means that they 
have a combined critical exponent that is equal to zero. Such amplitude ratios have been 
determined for a number of experimental systems and computed by using various the- 
oretical approaches like the e-expansion, perturbation theory in three dimensions fixed, 
high and low temperature series expansions and Monte Carlo simulations. A summary of 
results is given in refs. 0, 0]- Here we study universal amplitude ratios in the universality 
class of the three-dimensional Ising model with short range interactions, which is charac- 
terized by the Z2 symmetry of the order parameter. This universality class is supposed to 
be realized in a huge range of experimental systems: binary mixtures, uniaxial magnets 
or micellar systems; see 0, Isf. 

At finite values of the reduced temperature, power laws ([1]) are subject to corrections. 
For example the magnetic susceptibility behaves as 

X = C±\t\-'< {l + a±\tf + ht + ...) , (2) 

where 9 = uu = 0.524(4) [6]. The amplitudes C±, a± and b in general depend on the 
parameters of the system. Already in 1982 the authors of ref. [3] have demonstrated that 
for a model that interpolates between the Gaussian and the Ising model there is one value 
of the interpolation parameter, where a± vanishes. Renormalization group predicts that 
the zero of leading correction amplitudes is the same for all quantities. In the following we 
shall call a model with a± = an improved model. Studying improved models simplifies 
the accurate determination of amplitude ratios using Monte Carlo simulations or high and 
low temperature series expansions. Here we simulate the improved Blume-Capel model 
on the simple cubic lattice. For the definition of this model see the next section. Our 
main motivation to perform these simulations was to compute the energy density of the 
bulk system in a large range of inverse temperatures. This quantity is needed in our 
ongoing study of the thermal Casimir effect in the three-dimensional Ising universality 
class. Here we use the data generated for various quantities to update the estimates of a 
number of universal amplitude ratios. Computing universal amplitude ratios, we follow 
the strategy of js], [9|, where the spin-1/2 Ising model had been studied and more recently 
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id . where we had studied an improved model in the XY universahty class in three 
dimensions. 

Our results are essentially consistent with previous Monte Carlo studies js], 12, isf and 



the most recent analysis of high and low temperature series expansions IJ]. Typically we 
reduce the error bars by a factor of two to three compared with these studies. Estimates 
obtained by using field theoretic methods are typically by a factor of ten less precise than 
those obtained here. 

The outline of our paper is the following: First we define the model and the observ- 
ables that we have measured. Next we discuss the update algorithm and give details 
of our simulations. Using the data obtained, we extract numerical estimates for various 
universal amplitude ratios. These estimates are compared with those obtained in previous 
Monte Carlo simulations, from high and low temperature series expansions, field theoretic 
methods and experiments. Finally we conclude. 



II. MODEL AND OBSERVABLES 



The Blume-Capel model is characterized by the reduced Hamiltonian 



<xy> 



where the spin might assume the values G {—1,0,1}. The sites on the simple cubic 
lattice are denoted by x = (xq, xi, X2) with Xj = 0, 1, Lj — 1. In the following we shall 
consider lattices with L = Lq = Li = L2 in the high temperature phase and Lq = 2L, 
Li = L2 = L lattices in the low temperature phase. Throughout we consider periodic 
boundary conditions. The first sum in eq. ([3]) runs over all pairs of nearest neighbor 
sites < xy > on the lattice and /3 = l/ksT is the inverse temperature. The partition 
function is given hy Z = 'Yli{s} exp(— , where the sum runs over all spin configurations. 
In the following we shall consider a vanishing external field /i = 0. The parameter D 
controls the density of vacancies s^; = 0. In the limit D ^ —00 vacancies are completely 
suppressed and therefore the spin- 1/2 Ising model is recovered. \n d>2 dimensions the 
model undergoes a continuous phase transition for —00 < D < Dtri at a /3c that depends 
on D. For D > Dtri the model undergoes a first order phase transition. Refs. 1^ 17| give 



for the three-dimensional simple cubic lattice Dtri ~ 2.006,2.05 and Dtri = 2.0313(4), 
respectively. 

Numerically it has been shown that on the line of second order phase transitions there 
is a value D* of the parameter D, where leading corrections to scaling vanish. In ref. 
[ist we found D* = 0.641(8). One should note that little effort was made to estimate the 
systematical error due to subleading corrections to scaling. Recently we have determined 
D* = 0.656(20) [sl, where now systematical errors are taken into account. In \d\ we have 
simulated the model at D = 0.641 and D = 0.655 in the neighborhood of the critical 
point. Using a standard finite size scaling analysis we find 

/3c(0.641) = 0.38567122(5) (4) 
/3c(0.655) = 0.387721735(25) (5) 

as estimates of the inverse critical temperature. We also find that the amplitudes of 
leading corrections at D = 0.655 are reduced by at least a factor of 30 compared with the 
spin- 1/2 Ising model. 
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A. The energy density and the specific heat 



Here, we define the energy density as minus the derivative of the reduced free energy 
density with respect to ^ 



(6) 



\<xy> 



where V — L0L1L2. The specific heat is the derivative of the energy density with respect 
to /3. One finds 



C = — = - 
d/3 V 



\,<xy> 



S^Sy 



\<xy> 



(7) 



B. The magnetic susceptibility and the second moment correlation length in 
the high temperature phase 



The magnetic susceptibility x and the second moment correlation length ^2nd are de- 
fined as 



and 



where 



s2n(i 



X/F-I 



(9) 



(10) 



is the Fourier transform of the correlation function at the lowest non-zero momentum. 
In our simulations in the high temperature phase, we have measured F for the three 
directions A; = 0, 1, 2 and have averaged these three results. 



C. The magnetization, the magnetic susceptibility and the correlation length 
in the low temperature phase 

The magnetization in presence of a magnetic field is defined by 

m{KL) = ]^(Y.s}j , (11) 

where we assume, for simplicity, a fixed ratio Lq/L with L = Li = L2. The spontaneous 
magnetization is then defined as 

m(0, 00) = lim lim m{h, L) , (12) 
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where first the thermodynamic limit is taken. In a Monte Carlo simulation it is too 
cumbersome to follow this route. Note that m(0, L) at a finite value of L is however 
exactly zero for symmetry reasons. 

To avoid this problem, Binder and Ranch 1^ proposed the following definition: 



mRMs{0,L) = — 



\ 



(13) 



Here, following eqs. (20,21) of [20 



we use 



mABs{0,L) = — 



(14) 



which in the low temperature phase converges faster than muMsi^, L). 
The connected two-point correlation function is given by 

G{x,y) = {sxSy) - {sx){sy) . (15) 

In the low temperature phase, for /i = we replace eq. (ITBl) . using eq. (IT^ . by 

Giow{x, y) l^^o = {s^Sy) - m^ss(0, L) . (16) 

In order to project to zero-momentum states of the transfermatrix, we consider the 
correlation function 

G{r) = {SoSr)-{S^){S,) . (17) 

of time slices 



S. 



1 



Xo 



V LiL-2 



Xl,X2 



'^{xo,Xi,X2) 



Note that with this normalization, the correlation function has a finite thermodynamic 
limit as I/i, L2 — )■ 00. In the low temperature phase, for vanishing external field h = we 
replace {So){Sr) by LiL2m%s{Q, L). 

The magnetic susceptibility can be written as 



r=— 00 

The effective correlation length is given by 



(19) 



-1/ln 



G{r + 1) 
G(r) 



(20) 



The exponential correlation length is defined as C,exp = hnir-s-oo ^e//('")- Since the trans- 
fermatrix is positive and symmetric, ^e// approaches C,exp monotonically from below. The 
second moment correlation length is defined by 



S2nd 



1^2 

2dx 



(21) 



where d = 3 is the dimension of the system and 

oo 

fi2 = dY, r^G{r) . (22) 

r=— oo 

Note that in the thermodynamic hmit, the definitions f l9|2T]) become equivalent. In the 
low temperature phase we have computed x and fi2 by using eqs. fll9|22p . respectively, in 
the following way: Up to a certain distance R we have used G{r) computed directly from 
the configurations that we have generated. Since the relative statistical error increases 
exponentially with the distance r, for r > Rwe have used instead 

G(r) = G(i?)exp(^-^-A^ . (23) 

In the following analysis we have used the data obtained by choosing R ^ A^eff{R). We 
have checked that these results are consistent with those obtained for R ^ 3^e//(-R)- 



III. THE SIMULATIONS 

A. The Monte Carlo algorithm 



Analogous to [21|, we have simulated the Blume-Capel model using a hybrid of local 
updates and single cluster updates In the high temperature phase we have used as 
local update the heat-bath algorithm. With the local update we run through the lattice in 
typewriter fashion. Running through the lattice once is called one sweep in the following. 
After two heat-bath sweeps we perform a certain number Nd of single cluster updates. 
We have chosen to be roughly one third of the number of lattice sites V divided 
by the average size of a cluster. In the following we shall denote two heat-bath sweeps 
followed by N^i single cluster updates as one cycle of the update. In the high temperature 
phase we have used the cluster algorithm to compute improved estimators of the magnetic 
susceptibility and the second moment correlation length. 

In the low temperature phase we have used a local Metropolis update that is imple- 
mented in multispin coding technique [23]. Details of our implementation can be found 
in [6|. Here, after ten sweeps of the local update, we performed Nd single cluster up- 
dates. Also here we have chosen N^i to be roughly one third of the number of lattice 
sites V divided by the average size of a cluster. Here we denote ten sweeps of the local 
update followed by single cluster updates as one cycle. In the low temperature phase 
we did not use cluster improved estimators since they do not reduce the statistical error 
significantly in this phase. 

As random number generator we have used the SIMD-oriented Fast Mersenne Twister 
algorithm [2^ . 



B. Simulations in the high temperature phase 

First we have checked which lattice sizes are needed to keep the deviation from the 
thermodynamic limit smaller than the statistical error of the observables that we measure. 



6 



Based on finite size scaling theory [25|, we expect that the L-dependence of a singular 
observable A is given by 

A{L, (3) = A(oo, /3) X [1 + gA{L/m)] , (24) 

where we have ignored corrections to scaling. In the absence of a massless mode, as it is 
the case here, 

gA{L/O^CAexp{-L/0 (25) 

for large values of L/^. At D = 0.655 and (3 = 0.372 we have simulated lattices with the 
linear lattice sizes L = 18, 20, 22, 24, 26, 28, 32, 40 and 48 to check the size dependence of 
the observables. In total these simulations took about 6 days of CPU time on single core of 
a Quad-Core Opteron(tm) 2378 CPU (2.4 GHz). Our results are summarized in table HTl 
We give the number of update cycles (stat), the number of single cluster updates N^i per 
update cycle and the estimates of the energy density E, the magnetic susceptibility x and 
the second moment correlation length C,2nd- We have fitted these results with the ansatz 
A{L) = A{oo) + CAexp{—L/^exp), where we have taken ^^xp = 3.09394(13), which is the 
result for the exponential correlation length that we have obtained for L = 48. Skipping 
the data for L = 18 we get an acceptable x^/d-o.f. for all three quantities. The result 
for the thermodynamic limit y4(oo) is given in the last row of table [Til The correction 
amplitudes are ce = 0.114(3), = —21.2(1.1) and c^2nd ~ —1-28(8). This means that the 
deviation from the thermodynamic limit is of the same size as the statistical error that 
we have reached here for L/^ f» 10 for all three quantities that we have studied. 

In the following simulations we have chosen L ^ ^0^2nd throughout. Since for most 
of our simulations L is clearly larger than 10^2nd and the relative statistical error of the 
magnetic susceptibility and the second moment correlation length is larger than that of 
the results discussed above, we expect that deviations from the thermodynamic limit can 
be safely ignored. 

In the case of D = 0.655 we have simulated at 201 different values of /3 starting from 
(3 = 0.25 up to /3 = 0.3872. For each value of f3 we have performed 500000 update cycles. 
At /3 = 0.3872 we have simulated an L = 300 lattice and find ^2nd = 26.698(7). In the 
case of D = 0.641 we have only simulated at 12 different values of (3 from /3 = 0.3827, 
where 6nd = 8.8993(7), up to /3 = 0.3849, where ^2nd = 20.859(3). Here, for the two 
smallest values of (3 we performed about 3 x 10^ update cycles and about 1.5 x 10^ for 
the larger ones. The simulations in the high temperature phase at Z) = 0.655 and 0.641 
together took about one year of CPU time on a single core of a Quad-Core Opteron(tm) 
2378 CPU (2.4 GHz). 

C. Simulations in the low temperature phase 

Also here we have checked which lattice sizes are needed to keep deviations from the 
thermodynamic limit sufficiently small to be safely ignored. To this end, we performed 
simulations a.t D = 0.655 and (3 = 0.405 using the linear lattice sizes L = 12, 16, 20, 24, 
30 and 40. These simulations took about three weeks of CPU time on a single core of 
a Quad-Core Opteron(tm) 2378 CPU (2.4 GHz). Our results for the various observables 
are summarized in table [TTl Starting from L = 20 the results are consistent among each 
other. 
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TABLE I: The energy density E, the magnetic susceptibihty x ^nd the second moment correla- 
tion length ^2nd at D = 0.655 and (3 = 0.372 for various linear lattice sizes L. Furthermore we 
give the number of update cycles (stat) and the number of single cluster updates Nd per update 
cycle. In the last raw we give the result of our extrapolation to the thermodynamic limit, as 
discussed in the text. 



L stat/106 N,i E X 6nd 



18 


50 


200 0.492029(5) 


24.4851(18) 3.07593(13) 


20 


40 


300 0.491812(5) 


24.5198(17) 3.07794(13) 


22 


35 


400 0.491721(4) 


24.5405(16) 3.07910(12) 


24 


30.9 


500 0.491678(4) 


24.5469(15) 3.07945(11) 


26 


32.4 


700 0.491658(3) 


24.5502(13) 3.07963(10) 


28 


20 


1000 0.491643(4) 


24.5532(14) 3.07985(11) 


32 


10 


1000 0.491635(5) 


24.5538(18) 3.08000(13) 


40 


7.5 


2500 0.491633(4) 


24.5532(14) 3.07994(11) 


48 


5 


2000 0.491633(4) 


24.5545(15) 3.08008(12) 



0.4916314(17) 24.5549(7) 3.08000(5) 



TABLE IL Estimates for the energy density E, the magnetization m, the magnetic susceptibility 
X, the second moment correlation length ^2nd ^nd the exponential correlation length ^^xp at 
D = 0.655 and /? = 0.405 for various linear lattice sizes L. The results for x sind ^2nd s,re 
obtained with R = A and correspondingly, S,exp is approximated by ^^ffi!^), since A£,exp ~ 6. 
Furthermore we give the number of update cycles (stat) that we have performed. This number 
includes the factor 64 of copies of the system that we have simulated. In all cases we have 
performed two single cluster updates for each update cycle. 



L stat/10^ ^ m X §2nd ^exp 



12 


640 


1.063260(12) 0.514811(5) 4.4891(4) 1.5483(3) 1.6056(6) 


16 


320 


1.063467(11) 0.514997(5) 4.4317(4) 1.5271(3) 1.5778(7) 


20 


320 


1.063487(8) 0.515010(3) 4.4281(4) 1.5258(3) 1.5759(6) 


24 


190 


1.063474(8) 0.515005(3) 4.4283(4) 1.5259(4) 1.5761(8) 


30 


130 


1.063482(7) 0.515008(3) 4.4284(4) 1.5258(4) 1.5762(8) 


40 


32 


1.063474(9) 0.515004(4) 4.4277(7) 1.5243(7) 1.5729(15) 



Below we followed the recommendation of [8,] and have used lattices with L > 20^2nd- 
Given the observations made here, this is surely a safe choice. 

In the case of D = 0.655 we have simulated 64 different values of (3 starting from P = 
0.3884, where ^2nd = 11.687(45) up to /3 = 0.42, where 6nd = 1.0293(7). In addition we 
have simulated at 85 values of /3 up to /9 = 0.60, where we have only measured the energy 
density. In the case of D = 0.641 we have simulated at 12 different values of P starting 
from /3 = 0.3866, where ^2nd = 9.589(19), up to /3 = 0.3899, where ^2nd = 4.3749(34). In 
total these simulations in the low temperature phase took about 10 years of CPU time 
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on a single core of a Quad-Core Opteron(tm) 2378 CPU (2.4 GHz). 



IV. UNIVERSAL AMPLITUDE RATIOS 



First we have computed the ratio of the amphtudes of the magnetic susceptibihty in 
the high and the low temperature phase. To this end one could fit the data for the 
magnetic susceptibility with an ansatz such as eq. for the data in the high and the 
low temperature phase separately. Using the results for C+ and C_ obtained this way one 
could compute the ratio C+/C_. Instead, following ref. jsj we use a different strategy. 
The amplitude ratio can be defined as 

^ = 1-4^^. (26) 

Following this definition, we have first calculated the the ratio x{t)/xi^t) finite values 
of t = f3c — (3. To this end, we have computed for the values f3iow, where we have 
simulated in the low temperature phase corresponding values (3high = 2(3c — (3iow Here 
we made no effort to simulate exactly at these values of f3high- Instead we interpolate 
between the values that we have simulated. To this end, we take the /3_ < Phigh < P+ 
which are closest to f3high- We compute c(/3±) = x(/3±)(/3c — /3±)''', using 7 = 1.23719 
jsj. Then we linearly interpolate to get an estimate of c{Phigh)- Finally we compute 
xiPhigh) = c{(3high){(3c — (3high)~'^- Siuce we have simulated a large number of /3- values, 
the systematical error introduced by this interpolation should be negligible. In figure [T] 
we have plotted our results for t) as a function of t. From RG theory we expect 

that 



, af + ht + cf + dt'*' + ef + ... , (27) 
x(— t) C_ 

where a and in particular d should be small here, since D = 0.655 is a good approximation 
of D*. Following ref. j26[ 6' = 1.05(7). Therefore the term ct^' can be hardly discriminated 
from the analytic correction bt. The term ef^ is caused by the analytic background of the 
magnetic susceptibility. Figure [T] suggests that x(^)/x(~^) is essentially a linear function 
of t. Therefore, the terms explicitly given in eq. fl27|) should be sufficient to fit xi^) / x{~t) ■ 
In particular we have fitted our data with the ansatze 

^^^^ -^ + bt (28) 



xi-t) 

^^^^ - ^ + at' + bt (29) 



x{-t) c_ 

^^^^ - ^ + at' + bt + (30) 



x{-t) c. 

using 6 = 0.524 and 7 = 1.23719. Using these ansatze we have performed a large number 
of fits. Below we give the results of those fits that include a maximal number of data 
points under the condition that x^/d.o.f. is close to one. The fits are done using our 
data for D = 0.655 if not stated otherwise. Fitting the data that satisfy t < 0.005 
with the ansatz we get C+/C_ = 4.7089(14), b = 12.0(4) and xVd-o.f.= 24.6/20. 
Fitting all available data for D = 0.641 with ansatz we get = 4.7145(21), 
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■ 0.005 0.01 0.015 0.02 

t 

FIG. 1: We plot the ratio xi^) / xi~'t) & function of t = /3c — /5 computed from our data for 
D = 0.655. For a discussion see the text. 



b = 11.4(9) and xVd-o.f.= 9.6/8. Using ansatz ([29]), fitting data with t < 0.019 we get 
C+/C- = 4.718(2), a = -0.43(6), b = 15.6(3) and xVd.o.f.= 55.8/53. Fitting the data 
that satisfy t < 0.022 with the ansatz (jSO]) we get C+/C_ = 4.712(6), a = -0.05(30), 
b = 8.7(4.9), e = 12.3(8.2) and xVd.o.f.= 66.6/57. As our final result we quote 

^ = 4.713(7) , (31) 

which is chosen such that it covers all results, including their error bars, of the fits quoted 
above. We have estimated the error due to the uncertainty of Pc by redoing some of the 
fits using ratios computed with (3high = 2(/3c + error) — (3iow We find that it is clearly 
smaller than the error quoted above. 

Next we have computed the amplitude ratio 

= . (32) 

/2nd,- *\0 C,2nd{—'t) 

The calculation is analogous to that of the ratio C^/C^. Therefore we abstain from giving 
details and directly quote our final result 

^^^ = 1.939(5) . (33) 

/2nd,- 
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0.02 



FIG. 2: We plot the ratio u{t) = 2,x{t) /il^^{t)m'^ {t) function of — t = /3 — /3c computed from 
our data for D = 0.655. For a discussion see the text. 



Next we have computed the RG-invariant quantity 

Qc = ^-^^^ , (34) 

where B is the amphtude of the spontaneous magnetization in the low temperature phase. 
To this end, we have first evaluated r = ^2nd/ X for all values of (3 that we have simulated in 
the high temperature phase, where we have computed the statistical error of this combined 
quantity using the Jackknife method. Then we have computed 

= hm r(t)m'^(-t) (35) 
t\o 

in the same fashion as we have computed C^/C^ and f2nd,+/ f2nd- above. As our final 
estimate we quote 

Qc = 0.3293(2) . (36) 
The renormalized coupling in the low temperature phase is given by 

u = -pi --r = lim — — - — r— — . (37) 

fU-B' *^oe|„,(t)m2(t) ^ ^ 

In figure [2] we plot 3x(i)/^2nd('^)"^^('^) a function of —t. Since only quantities in the 
low temperature phase are involved, there should be no analytic correction. However, 
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since 6' ^ 1 we kept a term bt in our ansatze. Based on various fits we arrive at the final 
estimate 

u* = 14.08(5) . (38) 

Now let us consider the ratio C,exp/^2nd- It turns out that it is difficult to determine 
the exponential correlation length accurately in the low temperature have. The time slice 
correlation function behaves as 

G(r) = ci exp(-r/ei) + Ca exp{-r/^2) + ■■■ ■ (39) 



Since the ratio ^1/^2 = 1.83(3) j27| is rather small, the effective correlation length ^eff, 
eq. f l20|) . converges only rather slowly to ^exp- On the other hand, the relative statistical 
error of G{r) increases exponentially. Therefore very large distances that are needed 
to get a small deviation of ^e// from C,exp are not accessible. As compromise, we have 
taken C,eff{R) with R ^ 4^e//(-R) as our final estimate. To check the systematical error 
introduced this way, we have compared our result with that for R ^ 3^e//(-R). As our 
final estimate we quote 

.020(5) . (40) 



fexp,- 



2nd,- 



Here the error should cover both the systematical deviation of ^e// from S,exp as well as 
systematical errors due to subleading corrections that are not included in our fits. In order 
to get more precise results for ^^xp and as a consequence for f exp,-/ f 2nd,- a variational 



analysis of a large set of correlation functions, as it has been done in ref. [27| would be 
useful. Furthermore the method of [28| to reduce the variance of correlation functions 
could help to compute G{r) accurately at large distances r. 



A. Ratios that involve the specific heat 

In order to compute amplitude ratios that involve the specific heat, we have analyzed 
our data for the energy density which can be accurately determined in the simulation. 
In the case of the energy density we have to separate the analytic background and the 
singular part, which is needed here. In the neighborhood of the critical point, the energy 
density behaves as 

E = E, + E, , (41) 
where the analytic background can be Taylor expanded around the critical point: 

E,{P) = + CUP - Pc) + dUP - Pc? + ... ■ (42) 

The singular part is given by 

Es = a±\t\-'^ yi{l + h±\t\^ + ct + ...) . (43) 

In a first step we have analyzed data generated in relation with [Gj] for cubic systems 
with a linear size up to L = 360 and periodic boundary conditions directly at the critical 
point. At the critical point the energy density behaves as 

E = Ens + aL-'^^^'" [l + cL~'' + ...) (44) 
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and the specific heat as 

C = Cns + bL-^+^/''{l + dL-'^ + ...) . (45) 

Here we perform fits fixing u = 0.63002(10) as obtained in [(I]. Our final estimate is 
taken from fits without any correction term and all lattice sizes that are larger or equal 
to Lmin = 24 taken into account. Systematic errors are estimated by performing fits that 
include corrections with an exponent that is either 0.832, 1.6 or 2. For D = 0.655 we get 

Ens = 0.602111(1) + 0.006 X (z/ - 0.63002) + 42 x (/3, - 0.387721735) (46) 

and 

a = 1.7490(5) + 14 x (z/ - 0.63002) - 1800 x {(3^ - 0.387721735) . (47) 

We have redone the fits with slightly shifted values of the input parameters u and (3c 
to obtain the dependence of Ens and a on these parameters. For the specific heat at 
D = 0.655 we get 

Cns = -19.1(1) - 1700 X (z/ - 0.63002) - 1300000 x (/J^ - 0.387721735) (48) 

and 

b = 25.30(5) + 1350 x [u - 0.63002) + 620000 x {(3^ - 0.387721735) (49) 

from fits with L^m = 64. The error is dominated by systematical errors that we have 
estimated from fits that include corrections to scaling. 
For the energy density at D = 0.641 we get 

Ens = 0.604870(2) + 0.01 x (z/ - 0.63002) + 41 x {13c - 0.38567122) (50) 

and 

a = 1.749(1) + 14 X (z/ - 0.63002) - 1800 x {/3c - 0.38567122) . (51) 
For the specific heat at -D = 0.641 we get 

Cns = -19.1(2) - 1700 X (z/ - 0.63002) - 1000000 x (/3, - 0.38567122) (52) 

and 

b = 25.3(1) + 1350 X (z/ - 0.63002) + 500000 x {Pc - 0.38567122) . (53) 

Next we have analyzed our data for the thermodynamic limit in the neighborhood of 
the critical point using the ansatz 

E{P) = Ens + Cns{f3 - f3c) + a±|/3 - Pct'^ + dns{f3 - f3cf + fe±|/3 - , (54) 

where Ens, Cns obtained above and ^c = 0.387721735(25) and a = 0.10994(30) [I are 
input parameters, while a-t, dns and b± are the 5 free parameters of the fit. Using the 
results of these fits we have computed A^/A^ = —a^/a_ and P = (1 — Aj^/ A_) / a, which 
depends less on the input value for a than A^/A_. 

Fitting all data for D = 0.655 in the interval [Pc - 0.004, (3c + 0.004] we get A+/A_ = 
0.53611(7), P = 4.2195(6) and x^/d.o.f. = 55.5/51 using the central values of the input 
parameters. For the interval [(3c - 0.0075, /3c + 0.0075] we get A+/A_ = 0.53614(3), 
P = 4.2192(3) and x^/d.o.f. = 165.7/93. As a check we have also fitted with an ansatz. 
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where we have added a term oc (/3 — f3c)^ compared with the ansatz The results for 
A^/A_ and P change httle compared with those given above. It turns out that the error 
of and P is actually dominated by the error induced by the uncertainty of our 

input parameters, Ens, Cns, Pc and a. In order to estimate this error, we have repeated 
the fits using shifted values of these input parameters. For example, we have replaced Ens 
by {Ens + error). 

In order to check for the effect of leading corrections to scaling, we have fitted all our 
data at D = 0.641 using the ansatz dMD- We find A+/A_ = 0.53624(11), P = 4.2183(10) 
and x^/d.o.f. = 22.3/19 using the central values of the input parameters. This means 
that the results obtained ai D = 0.641 and D = 0.655 are fully consistent. 

We arrive at the final estimates 

^ = 0.536(2) , (55) 

where the error is dominated by the uncertainty of a, followed by the uncertainty of Cns- 
In contrast 

P = 4.22(1) (56) 

depends much less on the value of a. Its error is dominated by the uncertainty of Cns- 
This different behavior of Aj^/A^ and P is actually much more important in the case of 
the XY-universality class, where a is close to zero and therefore the relative accuracy of 
a is much smaller than in the present case. 
In order to compute the quantities 

Qg,+ = (57) 

and 

Q^,_ = (58) 
we have approximated the singular part of the energy density by 

E,it) = Eit) - Ens - Cnst - (59) 

Then we have computed 

m = te2ncims{t) - (60) 

The quantity Q± is then given by 

Q+ = a(l -a)limg(t) (61) 
Q_ = a(l -a)limg(t) . (62) 

We have fitted our data with the ansatz 

q(t) = q* + at . (63) 

In the high temperature phase we find by fitting all data with t < 0.006 the result 
q* = 0.19412(3). We have redone this analysis with shifted values of Cns and Ens to 
estimate the effect on our result for q*. It turns out that the error is dominated by the 
errors induced by the uncertainty of Cns and Ens- We have also redone the analysis using 
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TABLE III: Results for universal amplitude ratios obtained by high and low temperature series 
expansions of three different improved lattice models and Monte Carlo simulations of the 
spin-1/2 Ising model and the improved 4>^ model fl^]. In all these cases a simple cubic 

lattice has been studied. 



Ref. 


A+/A. 


C+/C- 


Q+ Q- 

J2nd,~ 


Qc 


here 


0.536(2) 


4.713(7) 1.939(5) 


1.020(5) 0.01899(10) 0.00487(2) 14.08(5) 


0.3293(2) 


[30] 






14.25(12) 




M 


0.532(3) 


4.76(2) 1.956(7) 


0.01880(8) 0.00472(5) 


0.3315(10) 


m 




4.75(3) 1.95(2) 


1.017(7) 14.3(1) 


0.328(5) 


[9] 


0.560(10) 








m 




4.756(28) 1.935(14) 




0.326(3) 


[13J 


0.532(7) 









our data for D = 0.641. We get an estimate for q* that is fully consistent with that for 
D = 0.655. We arrive at the final result 



Q+ = a(l - a)q* = 0.01899(10) . (64) 
Performing a similar analysis we arrive at 

Q_ = 0.00487(2) . (65) 



V. COMPARISON WITH RESULTS GIVEN IN THE LITERATURE 



A. Monte Carlo simulations and high and low temperature series 



In table 
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TTTl we confront our results with those of previous Monte Carlo simulations 
with a comprehensive analysis of high and low temperature series l4| and 



the low temperature series estimate of u* given in [30[. In [1^ a parametric representation 
of the equation of state has been used to obtain results for the critical isotherm and the 
low temperature phase from high temperature series. For an exhaustive overview of the 
literature see table 11 of 

In gg we have simulated the spin-1/2 Ising model on simple cubic lattices of a linear 
size up to L = 120 and L = 128, respectively. Also the authors of [l3| have simulated the 
spin-1/2 Ising model on the simple cubic lattice. They have simulated at a large number 



of /5-values in both phases of the model on lattices of a size up to L = 128. In 12| the 0"^ 
model on the simple cubic lattice has been simulated at A = 1.1, which is the estimate of 
A* obtained in ref. [29^. The authors have simulated lattices up to the size L = 120. In 
addition to simulations at a vanishing external field h = 0, they have simulated h at 
the critical temperature. This allowed them to compute additional universal amplitudes 
ratios that we do not discuss here. 

Essentially our results confirm those of the previous work. Even in the worst case, the 
deviation between our result and that of the other works summarized in table IIIII is less 
than three times the combined error. 
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TABLE IV: Results for universal amplitude ratios obtained by using the e-expansion (e) and 
perturbation theory in three dimensions (3D) fixed. For the definition of the amplitude ratios 
and a discussion see the text. 



Ref. Method C+/C_ u* Q, 



m 


e 




4.8 1.91 




[32] 


e 


0.55 


4.8 




[36] 


e 






0.01966(17) 


[33] 


e 


0.44 


4.9 


0.0223 


m 


e 


0.524(10) 






[35] 


e 


0.527(37) 


4.73(16) 




[39] 


3D 


0.541(14) 


4.77(30) 


0.331(9) 


m 


3D 






0.01968(15) 


[40] 


3D 


0.540(11) 






[38] 


3D 




4.72(17) 2.013(28) 


14.2 


[35] 


3D 


0.537(19) 


4.79(10) 




m 


3D 






0.0203 



B. Field theoretic methods 



In table II Vl we have summarized results obtained from the e-expansion and perturbation 
theory in three dimensions fixed. Mostly we have taken these results from table 12 of ref. 
Note that in ref. |35|] the field theoretic methods have been used in connection with 



a parametric representation of the equation of state. In ref. [39| A+C+/i?^ = 0.0594(11) 
is given. By using the value of given e.g. by [s^l Qc can be computed. Here we only 
report those amplitude ratios that we have computed in this work. For a comprehensive 
list of amplitude ratios see table 12 of ref. Essentially the field theoretic results are 
consistent ours, albeit their accuracy is clearly lower than ours. The errors for (5+ given 



by (36|, |37| seems to be underestimated. 



C. Experiments 

Here we just mention the results of two experimental works to give the reader an idea 
of the accu racy that can be reached. Studying a mixture of succinonitrile and water the 



authors of |42( found A+/A_ = 0.536 ± 0.005 and Q+ = 0.0187 ± 0.0013. Studying the 
antiferromagnet FeF2, the authors of [43] found A^^^/A^ = 0.53 ± 0.01 and C+/C_ = 
4.6 ± 0.02. In particular for A^/A_ the accuracy of the experimental studies is close to 
ours. The results of both studies are consistent with ours, confirming universality. For a 
comprehensive summary of experimental results see refs. 0, lEf. 
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VI. SUMMARY AND CONCLUSIONS 



We have simulated the Blume-Capel model on the simple cubic lattice aX D = 0.641 and 
0.655 for a large number of inverse temperatures /? in a neighborhood of the critical point. 
These values of D are close to D* = 0.656(20), where the amplitudes of leading corrections 
to scaling vanish. We have simulated lattices up to 300'^ in the high temperature and 
500 X 250^ in the low temperature phase. Throughout we have chosen the linear size 
L of the lattice such that L ^ 10^2nd in the high temperature and L ^ 20^2nd in the 
low temperature phase to avoid significant deviations from the thermodynamic limit. 
In the high temperature phase at D = 0.655 we have reached the correlation length 
^2nd = 26.698(7). Using the data obtained in these simulations we have extracted precise 
numerical estimates for a number of universal amplitude ratios. We carefully estimated 
systematical errors caused by subleading corrections. 

In table Hm we have summarized our results and compare them with previous estimates 
obtained from Monte Carlo simulations or from high and low temperature series expan- 
sions of lattice models. Our results are essentially consistent with but more precise than 
previous estimates. The same holds for the comparison with field theoretic methods. Also 
the accuracy of experimental results given in the literature is lower than ours. 
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